Homework 3

library(GGally)
library(tidyr)
library(plotly)
library(plyr)
library(dplyr) 
library(MASS)
library(vcd)
library(tidyquant)

\(1.\) Parallel Coordinates

\((a)\) Draw a parallel coordinates plot of the data in “ManhattanCDResults.csv” in the data folder on CourseWorks. (Original data source and additional information about the data can be found here: https://cbcny.org/research/nyc-resident-feedback-survey-community-district-results). Your plot should have one line for each of the twelve Manhattan community districts in the dataset.

cd <- read.csv('../data/ManhattanCDResults.csv')
cd[3:14] <- lapply(cd[3:14], function(x) as.numeric(sub("%", "", x)))
cd <- gather(cd, key="district", value="value", cd1:cd12)
cd <- cd[c(3,1,4)]
cd <- spread(cd, key=Indicator, value=value)

ggparcoord(cd, 
           columns = c(2:46), 
           groupColumn = "district", 
           scale = "globalminmax", 
           alphaLines = .7) + 
  coord_flip() + 
  xlab("Indicators") + 
  ylab("Approval")

\((b)\) Do there appear to be greater differences across indicators or across community districts? (In other words, are Manhattan community districts more alike or more different in how their citizens express their satisfaction with city life?

1. There appears to be a greater differences across indicators than across community districts.
2. Manahattan community districts are more a alike in how they express their satisfaction with the city life.

\((c)\) Which indicators have wide distributions (great variety) in responses?

1. Neighborhood parks, Neightborhood as a place to live, Neighbrohood Playgrounds, Availability of cultural activities are the ones that stand out as having a wide distribution in responses.

\((d)\) Does there appear to be a correlation between districts and overall satisfaction? In order words, do some districts report high satisfaction on many indicators and some report low satisfaction on many indicators or are the results more mixed? (Hint: a different color for each community district helps identify these trends).

1. There does appear to be a correlation between districts and overall satisfaction. There are some indicators that have higher variety in satisfaction than others, but ultimately there is less results that are mixed.

\(2.\) Mosaic Plots

Using the “Death2015.txt” data from the previous assignment, create a mosaic plot to identify whether Age is associated with Place of Death. Include only the top four Age categories. Treat Age as the independent variable and Place of Death as the dependent variable. (Hint: the dependent variable should be the last cut and it should be horizontal.) The labeling should be clear enough to identify what’s what, that is, “good enough,” not perfect. Do the variables appear to be associated? Describe briefly.

• (Assuming the “top four Age categories”, means the top age categories meaning, 1, 1-4, 5-14, 15-24, that being the age from the top.)

1. Place of death and age do appear to have some association. Very little.
2. The likelihood of dying at the Descendent’s home increases with age.
3. The likelihood of dying at the Medical Faility - Dead on Arrival appears to be higher for ages 1 or younger and decreases with age.
4. The associations aren’t strong. They’re subtle, but some do exist.

death <- read.csv('../data/Death2015.txt', sep='\t')
ages = c("1", "1-4", "5-14", "15-24")
death <- subset(death, Ten.Year.Age.Groups.Code %in% ages)
death$Ten.Year.Age.Groups.Code <- factor(death$Ten.Year.Age.Groups.Code, levels=ages)

Really annoying (I couldn’t increase the size of the chart to see the labels.)

par(mar=c(7,8,7,7))
just_labels <- c("left", "center", "center", "center")
labels <- list(gp_labels = gpar(fontsize = 14, fontface = 3))
mosaic(Place.of.Death~Ten.Year.Age.Groups.Code, death, keep_aspect_ratio=FALSE, labeling = labeling_border(rot_labels = c(8,0,0,0), just_labels = just_labels)) 

\(3.\) Time Series

\((a)\) Use the tidyquant package to collect stock information on four stocks of your choosing. Create a line chart showing the closing price of the four stocks on the same graph, employing a different color for each stock.

from <- today() - years(1)
tickers <- c("AAPL", "MSFT", "AMZN", "GOOG")
stock.prices <-  tq_get(tickers, get="stock.prices", from = from)
close.prices <- spread(stock_prices, key=symbol, value=close)

ggplot(stock.prices, aes(date, close, color=symbol)) +
    geom_line()

\((b)\) Transform the data so each stock begins at 100 and replot. Do you learn anything new that wasn’t visible in part (a)?

1. Apologies for the excessive code. I couldn’t find a better way of normalizing to a base-index of 100 so I wrote my own function. It would be awesome if you could show me a better way of doing this.
2. We learn that the performance of GOOG and MSFT performed better than depicted in the line chart above.
3. Indexing the data helps to better analyze the changes of price over time.
4. Resource: https://www.dallasfed.org/research/basics/indexing.aspx

# Function to transform each stock/symbol to begin at 100.
normalize <- function(symbol){
  index.value <- close.prices[[symbol]][1]/100
  ans <- list()
  for(price in close.prices[[symbol]]){
    ans <- append(ans, (price/index.value))
  }
  df <- data.frame(symbol=unlist(ans))
  colnames(df) <- c(symbol)
  return(df)
}

close.prices$AMZN <- normalize('AMZN')$AMZN
close.prices$AAPL <- normalize('AAPL')$AAPL
close.prices$GOOG <- normalize('GOOG')$GOOG
close.prices$MSFT <- normalize('MSFT')$MSFT

ggplot(close.prices) + 
  geom_line(aes(date, AAPL, color="AAPL")) + 
  geom_line(aes(date, AMZN, color="AMZN")) +
  geom_line(aes(date, MSFT, color="MSFT")) +
  geom_line(aes(date, GOOG, color="GOOG")) +
  xlab("") +
  ylab("Close Price")

\(4.\) Missing Data

For this question, explore the New York State Feb 2017 snow accumulation dataset available in the data folder on CourseWorks: “NY-snowfall-201702.csv”. The original data source is here: https://www.ncdc.noaa.gov/snow-and-ice/daily-snow/

\((a)\) Show missing patterns graphically.

Set Up Snow data

snow.url <- "https://www.ncdc.noaa.gov/snow-and-ice/daily-snow/NY-snowfall-201802.csv"
main.snow <- data.frame(read.csv(snow.url, skip=1))
snow <- as.matrix(main.snow)
snow[snow == "T"] <- 0.01
snow[snow == "M"] <- NA
snow <- as.data.frame(snow)
rownames(snow) <- snow$Station.Name

Tidy Snow Data

tidysnow <- snow %>%
  rownames_to_column("id") %>%
  gather(key, value, -id) %>%
  mutate(missing = ifelse(is.na(value), "yes", "no"))

Preserve key order

preserve.key.order <- function(){
  ordered.keys <- list()
  for(i in unique(tidysnow$key)){
    ordered.keys <- append(ordered.keys, i)
  }
  return(ordered.keys)
}
tidysnow$key <- factor(tidysnow$key, levels=preserve.key.order())

Column Missing Pattern (geom_tile())

ggplot(tidysnow, aes(x = key, y = fct_rev(id), fill = missing)) +
  geom_tile(color = "white") +
  theme(axis.text.x = element_text(size=20,angle=90),
        axis.text.y = element_text(size=20),  
        axis.title.x = element_text(size=40),
        axis.title.y = element_text(size=40,angle=90)) +
  xlab("") + 
  ylab("")

\((b)\) Is the percent of missing values consistent across days of the month, or is there variety?

1. There are very few stations that are consistent with reporting values.
2. There’s a lot of missing values and nothing is really consistent accross the days of the month; the values vary a lot.
3. It would be interesting to sort by location and see if there’s a pattern in missing data.

\((c)\) Is the percent of missing values consistent across collection stations, or is there variety?

1. It’s hard to tell if there’s a pattern in missing data. According to Missing Data by County, plot below, we could see that NYC under reports their values.
2. However, there’s a lot of points clustered around NYC, and of course some of those are unreported.
3. There’s a lot of variety in missing values across all the stations.
4. The percent of missing values is inconsistent across collection stations and has lots of variety.

snow.url <- "https://www.ncdc.noaa.gov/snow-and-ice/daily-snow/NY-snowfall-201802.csv"
snow <- data.frame(read.csv(snow.url, skip=1))
snow$na.count <- apply(snow, 1, function(x) sum(x=='M'))
ggplot(snow, aes(x = Longitude, y = Latitude)) +
  geom_point(aes(colour = na.count), alpha=.7) + 
  scale_colour_gradient(low = "grey", high = "green") + 
  ggtitle("Missing Data by County")

\((d)\) Is the daily average snowfall correlated with the daily missing values percent? On the basis of these results, what is your assessment of the reliability of the data to capture true snowfall patterns? In other words, based on what you’ve discovered, do you think that the missing data is highly problematic, or not?

1. There does not seem to be a correlation between the average snowfall and missing values.
2. Graphing both na_totals and daily average by each day we see a very weak negative correlation.
3. Negative correlation makes sense, given that if average snow fall is high then the total_na should be low.
4. The data is not that reliable in capturing true snowfall patterns.
5. The missing data is highly problematic in capturing true snowfall patterns.

snow.url <- "https://www.ncdc.noaa.gov/snow-and-ice/daily-snow/NY-snowfall-201802.csv"
main.snow <- data.frame(read.csv(snow.url, skip=1))
snow <- as.matrix(main.snow)
snow[snow == "T"] <- 0.01
snow[snow == "M"] <- NA
snow <- as.data.frame(snow)

Create a Function to get na_total and average for each day of the month.

get.na_avg <- function(){
  dates <- c(
    "Feb.1","Feb.2","Feb.3","Feb.4","Feb.5","Feb.6","Feb.7","Feb.8","Feb.9","Feb.10",
    "Feb.11","Feb.12","Feb.13","Feb.14","Feb.15","Feb.16","Feb.17","Feb.18","Feb.19",
    "Feb.20","Feb.21","Feb.22","Feb.23","Feb.24","Feb.25","Feb.26","Feb.27","Feb.28"
    )
  na_total <- list()
  average <- list()
  date_list <- list()
  for(date in dates){
    na_total <- append(na_total, as.numeric(sum(is.na(snow[[date]]))))
    average <- append(average, mean(as.numeric(snow[[date]]), na.rm = TRUE))
    date_list <- append(date_list, date)
  }
  temp_df <- do.call(rbind.data.frame, Map('c', date_list, na_total, average))
  colnames(temp_df) <- c('date', 'na_total', 'average')
  # set up the data frame
  temp_df$date <- factor(temp_df$date, levels=dates)
  temp_df$na_total <- as.numeric(as.character(temp_df$na_total))
  temp_df$average <- as.numeric(as.character(temp_df$average))
  return(temp_df)
}
na_avg <- get.na_avg()

Get Correlation

cor(na_avg$na_total, na_avg$average)

[1] -0.4602537

Plot Scatter

ggplot(na_avg, aes(na_total, average)) +
  geom_point() +
  geom_smooth()

Plot Results over Time for fun

ggplot(na_avg) + 
  geom_line(aes(date, na_total, color="na_total", group = 1)) + 
  geom_line(aes(date, average, color="average", group = 1)) +
  ggtitle("NA Total and Average Snow Fall For February")

Normalize Data Then Plot Again

# Function to transform each na_total and average to begin at 10.
normalize2 <- function(col){
  index.value <- na_avg[[col]][1]/10
  ans <- list()
  for(i in na_avg[[col]]){
    ans <- append(ans, (i/index.value))
  }
  temp_df <- data.frame(col=unlist(ans))
  colnames(temp_df) <- c(col)
  return(temp_df)
}
na_avg$na_total <- normalize2('na_total')$na_total
na_avg$average <- normalize2('average')$average

Normalized Plot

ggplot(na_avg) + 
  geom_line(aes(date, na_total, color="na_total", group = 1)) + 
  geom_line(aes(date, average, color="average", group = 1)) +
  ggtitle("NA Total and Average Snow Fall For February (Normalized)")